Disentangling flow and nonflow correlations via Bayesian unfolding of the 
event-by-event distributions of harmonic coeflicients in ultrarelativistic heavy-ion 

collisions 
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The performance of the Bayesian unfolding method in extracting the event-by-event (EbyE) dis- 
tributions of harmonic flow coefficients «„ is investigated using a toy model simulation, as well as 
simulations based on the HIJING and AMPT models. The unfolding method is shown to recover 
the input 112 — V4 distributions for multiplicities similar to those observed in Pb-f Pb collisions at the 
LHC. The effects of the nonflow are evaluated using HIJING simulation with and without a flow 
afterburner. The probability distribution of Vn resulting only from nonflow in HIJING is nearly a 
Gaussian and can be largely suppressed in the data-driven unfolding method used by the ATLAS 
Collaboration. The residual nonflow effects have no appreciable impact on the V3 distributions, but 
are found to affect the tails of the V2 and 114 distributions; these effects manifest as a small simul- 
taneous change in the mean and standard deviation of the u„ distributions. For the AMPT model, 
which contains both flow fluctuations and nonflow effects, the reduced shape of the extracted u„ 
distributions is found to be independent of pt in the low px region, similar to what is observed in 
the ATLAS data. The prospect of using the EbyE distribution of the harmonic spectrum aided by 
the unfolding technique as a general tool to study azimuthal correlations in high energy collisions 
is also discussed. 

PACS numbers: 25.75. Dw 



I. INTRODUCTION 

In recent years, the measurement of harmonic flow co- 
efScients w„ has provided important insight into the hot 
and dense matter created in heavy ion collisions at the 
Rclativistic Heavy Ion Collider (RHIC) and the Large 
Hadron Collider (LHC). These coefficients are generally 
obtained from a Fourier expansion of the particle distri- 
bution in azimuthal angle </>: 
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where w„ and $„ represent the magnitude and phase 
(event plane or EP angle) of the n^'^-order anisotropy. 
The coefficients v„ stem from the hydrodynamic response 
to various shape components or eccentricities e„ of the 
created matter. In non-central heavy ion collisions, £2 is 
mainly associated with the "elliptic" shape of the nuclear 
overlap region. However, £2 in central collisions and the 
other Cn terms arise primarily from fluctuations of the 
nucleon positions in the overlap region [1]. Model cal- 
culations suggest that w„ scales nearly linearly with e„, 
for n < 3 [2]. The proportionality constant is found to 
be sensitive to the properties of the matter, such as the 
equation of state and shear viscosity [3, 4]. 

Measurements of w„ coefficients for n = 1-6 [5-9] have 
established the importance of the fluctuations in the ini- 
tial state, and provided quantitative constraints on the 
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transport properties of the created matter. Most of these 
measurements focus on the Vn values averaged over many 
events, which mainly reflect the hydrodynamic response 
of the created matter to the average collision geometry 
in the initial state. Recently the ATLAS Collaboration 
made the first measurement of the probability distribu- 
tion of charged hadron u„, p{vn) [10], as well as the cor- 
relations between two or three EP angles of different or- 
der [11] in Pb-hPb collisions at ^/sWn =2.76 TeV. These 
new types of flow measurements provide additional in- 
formation on the nature of the fluctuations in the initial 
geometry and the non-linear effects in the hydrodynamic 
evolution [12- 16]. 

The study of these new observables requires new ex- 
perimental techniques. ATLAS extracts the p{vn) using 
a data-driven unfolding method based on the event-by- 
event (EbyE) distribution of cither the single-particle 
or the pair A0 [10]. The key step in the measurement is 
to determine a response function which connects the ob- 
served signal v!^^ with the true u„ signal, p{vf^^\vn)- The 
v^^ distributions are calculated using charged particles 
in the pseudorapidity range \r]\ < 2.5 and the transverse 
momentum range px > 0.5 GeV, which are then unfolded 
via the response function to obtain the true Vn distribu- 
tions, using the Bayesian unfolding technique [17, 18]. 
The response function mainly corrects for the smear- 
ing due to the effect of finite per-event charged particle 
multiplicity, but it also reduces possible nonfiow effects 
from short-range correlations, such as resonance decays, 
Bose-Einstein correlation and jets [3] . The efficacy of the 
method has been demonstrated using many data-driven 
cross-checks, however, a detailed Monte-Carlo study can 
provide a direct check of the performance of the method. 



as well as a measure of the influence of nonflow effects. 

In this paper, the performance of Bayesian unfolding 
is studied using the HIJING [19] and AMPT [20] event 
generators. We show that the unfolding procedure in- 
deed accounts for the effects associated with finite mul- 
tiplicity and nonflow. The structure of the paper is as 
follows: Section II gives a brief review of the unfold- 
ing method developed for the data analysis. Section III 
demonstrates the basic performance of the unfolding pro- 
cedure via a toy model simulation that includes only sta- 
tistical smearing. The influence of nonflow effect is stud- 
ied in section IV by applying the unfolding procedure 
to HIJING events with and without a flow afterburner. 
Section V gives a study of the Vn distributions in AMPT 
events. Section VI gives a summary of the results. Sec- 
tion VII discusses the prospect of applying the method 
to the study of the azimuthal correlation in high-energy 
collisions. 



II. UNFOLDING PROCEDURE 

The ATLAS analysis method starts with a Fourier ex- 
pansion of the azimuthal distribution of charged particles 
in a given event: 
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where v°^^ is the magnitude of the observed per-particle 
flow vector: v°^^ = {v°^l,v°^y). In the limit of large mul- 
tiplicity and the absence of nonflow effects, it approaches 
the true flow signal: v'^^ —J- v„- Hence one needs to de- 
termine the response function p{v°^^\vn) or p{v°^^\vn), 
which can be used to unfold these smearing effects. 

The response function is determined via a two- 
subevents method (2SE): the charged particles are di- 
vided into two equal subevents with rj > and ij < 0. 
The smearing effects are estimated from the distribu- 
tion of the difference of the flow vectors between the two 
subevents, Psub {(Vn^^)" — (^n°'"')^)i for which the phys- 
ical flow signal cancels. This distribution is observed to 
be well described by a 2D Gaussian with identical widths, 
Jjgj, , in both dimensions and hence the response function 
can be obtained as: 



p(4°'i4) = 2Psub(2[(«r^)'^-(t)r^)'']) (3) 

„ -- ,5n = 6,,j2. (4) 



where the factor of two accounts for the fact that the 
smearing in v°^^ is a factor of two narrower than that for 
^^^bs^a_^^^bs-jh^ Consequently the ID response function 
can be obtained by integrating out the azimuthal angle: 
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where Iq is the modified Bessel function of the first 
kind. Eq. 5 is known as the Bessel- Gaussian (BG) func- 
tion [21, 22], which is a very good approximation of the 
response function. In ATLAS data analysis, the actual 
distribution, obtained by shifting the measured 2D dis- 
tribution Eq. 3 to Vn and integrating out the azimuthal 
angle, is used in the unfolding. 

The Bayesian unfolding procedure from [17] is used to 
obtain the Vn distribution. In this procedure, the true w„ 
(cause "c") distribution is obtained from the measured 
v°^^ (effect "e") distribution and the response function 
Aji = p{ej\ci), using an iterative procedure: 
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The unfolding matrix M'^ is determined from the re- 
sponse function and initial guess of the true distribution 
cP (referred to as prior) , and the process is iterated. The 
V2 and V3 distributions were shown to converge after a 
few iterations, while more iterations (> 16) were needed 
for Vi. One of the goals of this paper is to verify the 
convergence of the unfolding via simulations. 

The ATLAS results also show that the response func- 
tion (Eq. 5) not only accounts for the smearing effect, but 
also for the nonflow contributions. The argument goes as 
follows: most of the nonflow effects arise from localized 
clusters each containing a few particles, and the number 
of such clusters is proportional to the total multiplicity. 
Therefore, the nonflow contributions to v°^^ = \v°^^\ 
is of a similar nature as that for the response function 
which is obtained from {v°^'^)°' — {v°^^)''. Thus, most 
of the short-range nonflow effects can be removed by the 
unfolding procedure. This can be verified directly via 
a Monte-Carlo event generator such as HIJING [19] as 
discussed later. 



III. TOY MODEL SIMULATION 

The toy model simulation serves as a proof-of-principle 
check for the unfolding method. The multiplicity of each 
event is generated according to a Gaussian distribution 
with a mean of (N) — 1070 and a width of (Jn ~ 84. 
This roughly simulates the distribution of the number of 
reconstructed charged particles for px > 0.5 GeV and 
I??] < 2.5 in 20-25% Pb-hPb collisions at ^/s]^ =2.76 
TeV in ATLAS, assuming a 70% reconstruction effi- 
ciency [23]. The input w„ signal is obtained by a Bessel- 
Gaussian function that is parameterized to the ATLAS 
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FIG. 1: The ATLAS preliminary Vn data [10] and parameterization by Bessel-Gaussian functions. These functions are used as 
the input distribution for Vn in the toy model and HIJING simulations. 
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FIG. 2: The performance of the Bayesian unfolding for V2 (left), V3 (middle) and V4 (right) in toy simulation. 



data as shown in Figure 1. Particles in a given event are 
implemented to have the same w„ value, which is sampled 
from the BG functions. The particles arc then randomly 
divided into two subcvcnts with equal multiplicity. The 
correlations of the pcr-particle flow vectors between the 
two subevents are then used to calculate the response 
function. About 2.5 Million events were generated; this 
is similar to the statistics used in the ATLAS data anal- 
ysis for the 20-25% centrality bin. 



Figure 2 shows the performance of the Bayesian un- 
folding for different number of iteration (iVitcr)- The re- 
sponse function (thin solid curves) are significantly nar- 
rower than the v';^^ distribution for n — 2 and n = 3, 
leading to convergence to the input distribution after a 
few iterations (see ratios for iViter = 8). However the 
convergence for n = 4 is slower due to its much smaller 
Vn signal. This overall performance is similar to what is 
observed in the data analysis [10]. 



IV. HIJING SIMULATION 



Understanding the nonflow effects in any flow mea- 
surement is a topic of intense debate over the last several 
years [3, 24-28]. In this paper, the nonflow correlations 
refer to statistical smearing due to a finite number of 
observed particles, short-range correlations among few 
particles such as resonance decays, Bose-Einstein corre- 
lation and fragments of individual jets, as well as long- 
range correlations not related to flow such as Global mo- 
mentum conservation (GMC) and fragments from di-jets. 
Note that most flow methods are not sensitive to pure 
statistical smearing, but it is listed as part of nonflow 
since it appears explicitly in the EbyE distributions. Also 
note that the GMC effect is a global correlation that 
mainly affects the vi, and the correlations between two 
back-to-back jets involve a few particles but can be long- 
range in 77, affecting mainly the even harmonics. These 
nonflow correlations and their fluctuations are naturally 
included in the HIJING event generator and has been 
used frequently in previous investigations of nonflow ef- 
fects [25, 26, 29, 30]. 



Since HIJING does not have physical flow, any non- 
zero v'^^ values must be solely due to nonfiow. However, 
we can also impose a physical flow signal on the particles 
using a flow afterburner, such that the generated events 
contain both the nonflow and flow contributions that fluc- 
tuate event- by- event. The Bayesian unfolding procedure 
can be used to extract the u„ distributions for compari- 
son with the input. This allows a quantitative estimate 
of how the nonflow effects and their EbyE fluctuations 
affect the final extracted w„ distributions. 

In this study, we generated 0.5-3 million HIJING 
events each at four fixed impact parameters, 6 = fm, 
b = A fm, 6 = 8 fm and 6 = 11 fm for Pb-|-Pb collisions 
at 2.76 TcV. The event multiplicity roughly corresponds 
to 0-1%, 5-10%, 20-25% and 50-55% centrality ranges in 
the data. The experimentally measured u„ distributions 
in these four centrality ranges arc parameterized indi- 
vidually by the BG functions. These functions arc then 
used as the input w„ distributions to be imposed on the 
HIJING events. The unfolding procedure is repeated for 
all four impact parameters, but the results for 6 = 8 fm 
are often used for illustrative purposes. In order to gain 
maximum sensitivity to the nonflow effects, all particles, 
including neutral particles, above 0.5 GcV and |?7| < 2.5 
arc used. 




A. Default HIJING (without flow) 

The goal of HIJING simulation without a flow after- 
burner is to understand the basic properties of the re- 
sponse function and the contributions of non-flow effects. 
Figure 3 shows the y-axis projections (the x-axis projec- 
tions are identical) of the difference of the per-particle 
flow vector between the two symmetric subevents (77 > 
and ?7 < 0), psub {{v,f''T " (4°''')^)- They are shown for 
n = 2-7 in events with 6 = 8 fm. These distributions are 
well described by a Gaussian, implying that the smearing 
effects and those due to nonflow short-range correlations 
are largely statistical. Figure 4 shows similar distribu- 
tions for various impact parameters for n — 2. Signiflcant 
non-Gaussian behavior is observed for peripheral colli- 
sions (i.e. 6 = 11 fm). In this case, the distribution can 
be described by the Student's t-distribution [10], which 
is a general probability distribution function for the dif- 
ference between two estimates of the mean from indepen- 
dent samples. The t-distribution approaches a Gaussian 
distribution when the event multiplicity is large. The dis- 
tributions shown in Figures 3 and 4 are used to obtain 
the response functions according to Eq. 3. 

Since Vn = 0, in the limit of pure statistical smearing, 
the response function reduces to the v°^^ distribution, 
and is a 2D Gaussian centered around the origin: 

P(C^^) = p(-ri^n = 0) = J-e-(^"^lV(2.S) (7) 
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The X or y width of the 2D response function or i5„ can 
be calculated from the central limit theorem as (see for 
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FIG. 3: The y-axis projection of the distribution of the 
difference between the per-particle flow vectors of the two 
subevents for HIJING events at b = 8 fm (i.e. projection 
of Psub ((Wn"""")" - («,?''")'')). Each panel shows the distribu- 
tion for one harmonic number. Lines indicates a fit witli a 
Gaussian function. 



FIG. 4: The y-axis projection of the distribution of the 
difference between the per-particle flow vectors of the 
two subevents for HIJING events (i.e. projection of 
Psub ((Wn"*"")" - («,?''")'')). Each panel shows the distribution 
at one impact parameter. Lines indicates a fit with a Gaus- 
sian function. 
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FIG. 5: Top panels: distributions of v,^^^ (diamond), response function obtained via standard 2SE method (squares) and 
response function obtained via two subevents with randomization of tracks in (p (circles). Bottom panels: The ratios of two 
response functions to the v^ " distribution. Each column shows the result for one harmonic number (n = 2-6 from left to 
right). 



example Ref. [21, 24]): 
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where N is total number of particles used in the event 
and the weight for each particle w is taken to be one 
(w = 1). This width is the limit expected for pure sta- 
tistical smearing. Additional nonflow correlations tend 
to increase this width. For example, if the N particles 
come from Nr resonance decays each with M number of 
particles with the same 4> angle, the distribution is still 
a Gaussian but with a larger width: 
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In the data analysis, the response function is obtained 
by 2SE method via Eq. 9, which has a slightly broader 
width: 



1 



2Na 



1 



1 



2A^(l-a2) 



J„(l + 0.5a^ 



(11) 



where a = {Na — Ni,)/{Na + Ni,) accounts for the fact 
that the two subevents may not have the same multi- 
plicity event-by-event. If the fluctuations between the 
two subevents are assumed to be random, one expect 
a^ w 1/N. The average a^ values are calculated di- 
rectly in HIJING for the four impact parameters, and 
they are found to be (a^) =0.05%, 0.07%, 0.14% and 
0.4%, respectively. This implies that the relative multi- 
plicity fluctuation between the subevents has negligible 
effects as long as the input Vn distribution and total mul- 
tiplicity N = Na + Nh remain unchanged. 



Eq. 9 is valid for events with fixed multiplicity. The 
spread of the particle multiplicity in the event class can 
further increase the width: 
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where aN is the standard deviation of the particle mul- 
tiplicity distribution. 

The argument in Eq. 10 can be generalized to other 
short-range correlations as long as the number of sources 
are proportional to the total multiplicity. On the other 
hand, long-range correlations such as di-jets can correlate 
between the two subevents (with each subevent contain- 
ing one jet), and hence their influence to the response 
function may differ from Eq. 10. Figure 3 shows that the 
distributions of (4°''")° - {v°^''f from HIJING arc con- 
sistent with a Gaussian, but their widths decrease with 
harmonic number n. The width for n = 7 is very close 
to the pure statistical fluctuation limit (see Figure 6). 
This behavior suggests that the nonflow contributions 
are largely statistical in nature, and their influence to w„ 
decrease with the harmonic number. 

Figure 5 compares the distribution of v°^^ and the re- 
sponse function obtained using the standard 2SE method 
(denoted as RF2se), as well as the response function 
obtained by randomly dividing the particles into equal 
halves (denoted as RFran)- This randomization proce- 
dure destroys all of the nonflow azimuthal correlations 
and thereby has a width consistent with Eq. 12. 

If RF2SE properly accounts for all nonflow effects, it 
should coincide with the v°^^ distribution. Figure 5 
shows that this is true for most harmonics, except for 
n ~ 2 and to a lesser extent for n = 4. These differences 
could be due to residual effects of di-jets, which affect 
mainly the even harmonics. It is also interesting to point 
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FIG. 6: The widths of the distributions vs n for v°^\, RF2SE 
and RFran in four centrality classes in HIJING without fiow. 
The solid lines indicate the values calculated via Eq. 9 using 
the average multiplicity, while the dashed lines indicates the 
values taking into account the multiplicity fluctuation Eq. 12. 
These distributions can be regarded as the power spectrum 
for nonflow in HIJING. 
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out that the three distributions approach each other for 
large n and reach the statistical limit given by Eq. 12 for 
n > 6. This behavior is quite natural, as n^'^-order har- 
monics are sensitive to correlations at an angular scale of 
27r/n. The harmonics at very large n eventually resolve 
the (j) angle of all particles individually, and hence depend 
only on the total multiplicity in the event. 

We found that all distributions in Figure 5 are well 
described by the radial projection of a 2D Gaussian cen- 
tered around the origin, i.e. Eq. 8. Hence they can be 
represented by a single parameter (5„. Figure 6 sum- 
marizes the 6n versus n for the three distributions for 
four impact parameters. The vifidths for RF2SE £^i'e found 
to agree with that for the v°^^ distributions, except for 
n = 2. The small difference seen for n = 4 may be associ- 
ated with di-jets that may have a non-zero V4 component. 
However, all three types of widths are found to approach 
the statistical limit at large n given by Eq. 12 (dotted 
lines). 

To quantify the nonflow effects that lead to the dif- 
ference between the response functions and the v°^'^, we 
unfold the v°^^ distribution using RF2SE and RFran- 
The unfolded distributions for n = 2 are shown in Fig- 
ure 7 for two centrality classes (6 := 8 fm and 6 = 11 fm). 
These distributions can be taken as the EbyE distribu- 
tion for the remaining nonflow: the distribution obtained 
with RFran represents the distribution of nonflow be- 
yond pure statistical limit, while the result obtained with 
RF2SE represents the distribution of nonflow that are 
correlated between the subevents (hence not included in 
RF2SE and remains after unfolding). The latter is nar- 
rower than the former as expected. Both types of non- 



FIG. 7: The distributions of V2'^'' (open diamonds) and the 
two response functions RF2SE (open squares) and RFran 
(open circles), together with the estimated residual nonflow 
distributions when the two response functions are used (fill 
squares and fill circles). The top panel and bottom panels 
shows the results for b = 8 fm and fc = 11 fm respectively. 
The curves represent a fit to Eq. 8. 
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FIG. 8: Per-particle nonflow estimated for HIJING when the 
response function RF2SE is used. The result are obtained 
by assuming quadrature difference (circles) and linear differ- 
ence (squares) between the Sn values for v^^'^ distribution and 

RF2SE. 



flow distributions arc nearly Gaussian, except in the tails. 



Similar unfolding exercises are repeated for higher-order 
harmonics. Nonflow signals are found to be quite small 
for n > 2, and their EbyE distributions could not be 
extracted reliably. 

The results shown in Figures 6 and 7 allow us to es- 
timate the contributions of nonflow relative to the flow 
signal, as shown in Figure 8. They are estimated either 
as the quadrature difference or as linear difference be- 
tween the 5n values for the v°^'^ distribution and RF2se- 
The former estimate is appropriate if the nonflow ef- 
fects contribute as independent smearing on top of the 
l/\/2N smearing. The latter estimate is more appro- 
priate if the nonflow is a overall shift of the v°^^^ dis- 
tribution. The reality probably sits between these two 
extremes. The worst case example is provided by values 
for 6=11 fm (the right most point in Figure 8): assum- 
ing the real flow signal is 0.1, the nonflow contribution is 
(0.018/0.1)2 =3.2% (quadrature) and 0.004/0.1=4% (lin- 
ear) . In general, the quadrature (linear) estimation gives 
a larger nonflow when the flow signal is small (large). 



B. HIJING with flow after-burner 

To understand how the nonflow effects influence the 
unfolding performance for extracting the real flow signal, 
HIJING events are modulated with a px dependent flow 
signal as: 
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where p{vn) denote the BG function obtained by fit- 
ting the ATLAS data as in Figure 1, w„(pt) is the pt- 
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FIG. 9: The response function using the standard subevents, 
RF2SE, for HIJING simulations without and with flow for 
n=2-4. The bottom part of each panel shows the ratio be- 
tween the two cases. 



differential f„ measured by the event-plane method [9], 
and Wn(po) is the integral v„ for px > 0.5 GcV from the 
event-plane method. The flow signal is implemented by 
rotating the (/) angle of the particles according to the pro- 
cedure discussed in Ref. [31]. The phases of «„, $„, are 
randomly distributed, thus there are no correlations be- 
tween different n. The amount of rotation is a smooth 
function of — ^„[31], and hence largely preserves the 
original non-flow correlations. The response functions 
RF2SE a-re unchanged by the implementation of the flow 
signal as shown in Figure 9. 

The unfolding performance for HIJING with flow is 
shown in Figure 10 for & = 8 fm. The unfolding proce- 
dure is found to always converge to a stable distribution 
after a few iterations. But this distribution can be dif- 
ferent from the input, as shown in the bottom panels of 
Figure 10. The agreement for RF2SE is systematically 
better than the result obtained with RFran since some 
nonflow effects are included in RF2se- The remaining 
differences from the input distribution are presumably 
due to residual nonflow effects not captured by RF2se- 
These nonflow effects lead to a broadening of the f 2 and 
W4 distributions, leading to sizable deviations in the tails 
of the distributions as shown by the ratios in the flgure. 
However, these residual nonflow effects have negligible 
effects on the f 3 distribution. 

Identical studies are also carried out for other impact 
parameters. The results for 5 = fm and 11 fm are 
shown in Figure 11 and Figure 12, respectively. The per- 
formance for h = 11 fm is noticeably poorer, i.e. the final 
distributions differ significantly from the input w„ distri- 
butions in the tails. This is a natural consequence of the 
influence of the residual nonflow effects not captured in 
the response functions. 

Figure 13 compares the mean value (un) and the 
standard deviation a^^ between the input distribution 
and the unfolded distributions obtained with RF2SE £^nd 
RFran- The values from the unfolded distributions are 
alway larger than those for the input distributions, sug- 
gesting that the effects of nonflow always increase the 
measured («„) and a^^. The v^ results for RF2SE agrees 
with the input within 4% for all four impact parameters; 
the agreement for U3 is better than 1% for 6 = 0, 4 and 
8 fm, and better than 3% for & = 11 fm; the agreement 
for W4 is better than 2% for 6 = and 4 fm. and worsen to 
6% for h — 11 fm (consistent with the simple estimation 
in Figure 8). 

Most of the deviations from the input values are corre- 
lated between (w„) and (t„„ , which iinplics that the main 
difference between the input w„ distribution and the un- 
folded Vn distribution in Figures 10-12 can be absorbed 
into a simultaneous change of {v^^ and cr„^. The shape 
of the unfolded distribution and the shape of the input 
distribution, when rescaled to have the same (wn), are in 
fact found to be quite similar. 

Thus far, our discussion has focused mainly on n > 2. 
In principle, the analysis procedure can also be ap- 
plied to measure the rapidity-even dipolar flow v\ , which 
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is found to be large and dominates over the rapidity- 
odd directed flow signal at LHC energies [9]. How- 
ever the vi analysis is sensitive to global momentum 
conservation effects, and a p^ and 77 dependent weight, 
w{pT,il) =PT- {PT)i'n)/{PT){v) [32, 33], needs to be ap- 
plied to each particle while calculating the response func- 



tions. We defer this topic to a future study. 
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AMPT MODEL SIMULATION 



A Muti-Phase Transport model (AMPT) [20] has been 
used frequently to study the higher-order u„ associated 
with e„ in the initial geometry [34-36]. It combines the 
initial fluctuating geometry from HIJING and final state 
interaction via a parton and hadron transport model. 
The parton and hadron transport is responsible for trans- 
forming the e„ into the momentum anisotropy. The 
AMPT model naturally contains various nonflow effects 
present in HIJING, but these nonflow effects can also mix 
with the collective flow during the transport process. In 
many aspects, the AMPT model provides a more realistic 
framework for understanding the role of the nonflow and 
the hydrodynamic response of the matter to the initial 



geometry. 

The AMPT sample used in this study was generated 
for 6 = 8 fm AuH-Au collisions at RHIC energy of 
y/SNN = 0.2 TeV, containing about 0.8 million events ^. 
The AMPT events were generated using the parameters 
that were tuned to reproduce reasonably the experimen- 
tal pt spectra and Vn data at RHIC. They were gener- 
ated with the string-melting mode with a total partonic 
cross-section of 1.5 mb and strong coupling constant of 
as = 0.33 [35]. 



^ Due to the large amount of resources required for the AMPT 
simulation, we don't have enough statistics of Pb-f-Pb collisions 
at the LHC energy for this study. 



11 



,0° 



♦ 



° V^=0.2TeV 

h0.5<p^<1 GeV 
hP^>1 GeV 
^p >0.5GeV 



0» 



0.1 0.2 0.3 

Vp 



10= 
10' 


/ \ 






10' 


r =■♦ - 




, , , , 1 , , . , 1 , , , , 1 . , , , 1 , 



0= 


■ '''!^2i*6to, 


1 


0' 


r ♦ n 


_ 




^ n 






° ♦♦ °o 


: 




^ 


- 


n^ 


•♦ p 


_ 




^0.5<p^<1GeV 0? 


: 




-e-P^>1 GeV ' 


- 




-»-p^>0.5GeV ij, ijil 


r 



0.1 0.2 0.3 0.4 

V, rescaled 



c 




0.1 0.2 

V3 


0.3 


0= 


/" 


X 






:r 

n 


\ 


■ 


0' 
0^ 


i 




: 






i 


1 1 1 




0.1 0.2 

V, rescaled 






0.05 


O.'l 


0.15 










V4 






lO*- 


■/ 


-\ 






: 


10' 


f 

1 


> 


\ 




1 


Itf 


- 


1 1 1 1 1 


\ 
\ 

1 1 1 1 . 


1 1 1 1 


- 



0.05 0.1 0.15 

V, rescaled 



FIG. 16: Top panels: The unfolded distributions for i)„ in & = 8 fm AMPT events at RHIC energy for particles in pt > 0.5 
GeV, 1 > Pt > 0.5 GeV and px > 1 GeV ranges. Bottom panels: same distributions but rescaled horizontally so the (vn) 
values match that for pt > 0.5 GeV range. 



Figure 14 shows the distributions of v°^^ and the two 
types of response functions (RF2SE and RFran) for n ~ 
2 — 6. The difFcrcnccs between the two response functions 
and their dependence on n are quite similar to what is 
observed in the HIJING simulation (Figure 6). The large 
difference in dn between v°^^ and the response functions 
are due to significant genuine collective Vn signal up to 
n = 6. This is clearly revealed in the bottom right panel 
of Figure 14. In this case, the (5„ for v°^^ is obtained via 
fitting to the BG function: 



p{x) 



6^ 



RP 

(52 



(14) 



v^^ is the anisotropy associated with the average shape 
of the overlap region, while (5„ has contributions from 
both nonflow and geometry fluctuations, u^^ is allowed 
as a free parameter for n = 2 but is fixed to be for 
n>2. 

The unfolding performance are shown in Figure 15, the 
influence of nonflow is reflected by the difference between 
the unfolded results for the two response functions. The 
difference is small for V2 , but become quite sizable in the 
tails for W3 and V4. 

One important finding in the ATLAS analysis is that 
the shape of the Vn distributions is nearly independent of 
the Pt of the particles used (in the low px region). This 
is checked directly in the AMPT simulation. The results 
in Figure 16 suggest an approximate scaling, although 
significant deviations arc observed for V2 in the tails. 



VI. SUMMARY 

The performance of the Bayesian unfolding method, 
used by the ATLAS Collaboration to extract the w„ prob- 
ability distributions for n = 2 — 4, is evaluated with a toy 
model simulation, as well as with HIJING and AMPT 
Monte-Carlo model simulations. In the absence of non- 
flow effects, the unfolding method is demonstrated to 
converge reliably to the input distributions for V2 — f 4 for 
the typical event multiplicity used in the ATLAS anal- 
ysis. The effects of the nonfiow are evaluated using HI- 
JING simulation with and without a flow afterburner. 
The nonflow effects in HIJING are found to be largely 
statistical in nature, which is consistent with them being 
associated with independent emitters each containing a 
flnite number of particles. The power spectrum for non- 
flow are found to decrease rapidly with n and approach a 
purely statistical limit given by ^J {1/2N) where N is the 
number of particles used in the event. The majority of 
these nonflow effects are included in the data-driven re- 
sponse function and hence are removed in the unfolding 
procedure. The probability distribution of the residual 
nonflow, attributed partially to the away-side jets, is ex- 
tracted for n = 2 and found to be largely Gaussian. The 
influence of these residual nonflow effects arc found to 
be significant in the tails of the V2 and W4 distributions 
and for the peripheral collisions. However this influence 
can be absorbed in to a simultaneous change, by a few 
percent, of the mean and the width of the v^ distribu- 
tions, such that the reduced shape is largely invariant. 
Thus the data-driven unfolding procedure seems to be 
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rather robust in removing the nonflow contributions and 
recovering the true long-range flow correlations. 

The unfolding technique is also applied to AMPT 
transport model simulation which contains both the flow 
fluctuations and nonflow. The Vn distributions are ex- 
tracted in two pt ranges in the low px region. When 
these distributions are rescalcd to the same mean values, 
the adjusted shapes are found to be nearly the same, very 
similar to the observation of the ATLAS measurement. 
This flnding is also similar to recent event-by-event hy- 
drodynamics model calculations [15, 16]. 



VII. EBYE DISTRIBUTION OF HARMONIC 

COEFFICIENT AS A TOOL FOR STUDYING 

AZIMUTHAL CORRELATIONS IN 

HIGH-ENERGY COLLISIONS 

Although we discussed the calculation of EbyE har- 
monic coefficients v'^^ and the Bayesian unfolding in the 
context of studying flow and nonflow in heavy ion col- 
lisions. This method can also be used as a general tool 
to understand the azimuthal correlations in other high 
energy experiments such as proton-proton or proton- 
nucleus collisions, as well as for events selected with par- 
ticular features, such as events with multiple jets or large 
missing Ex- Such studies are still in progress, neverthe- 
less but we briefly describe the idea as follows. 

Each collision produce a distribution of particles in az- 
imuth with a per-particle vf^^ spectrum calculated via 
Eq. 1. However this spectrum in a single event can not be 
used directly as it is dominated by the statistical "noise" 
due to the finite number of particles. This problem can be 
circumvented by calculating the distribution of vf^^ from 
many events. Each event class of interest can have its 



own set of v'^^ distributions, one for each n. The statis- 
tical noise for these distributions can then be estimated 
from the large n limit (in practice, n = 10 typically is 
enough, see Figure 6): 



Poo(«) = lim pivT). 



(15) 



This distribution can also be obtained by randomizing 
the (j) angle of all tracks (hence destroy all azimuthal cor- 
relations) . This distribution gives a baseline distribution 
in the absence of any azimuthal correlation. The distri- 
bution of lower-order harmonics is a convolution of this 
"noise" distribution with a probability distribution for 
genuine correlation p(s„): 



p{vT)=Poo{v)®p{Sn). 



(16) 



where the s„ is a per-particle measure of the correlation 
strength at angular scale 2-k /n. Standard unfolding tech- 
nique can then be used to solve for p{sn) using the Poa as 
the response function. The distributions of p(s„) in prin- 
ciple should be sensitive to the nature of the azimuthal 
correlations of the underlying physics processes in each 
event class. 

If one also wants to distinguish correlations that are 
long range in rj from those that arc short-range, one 
should use the response function obtained from the two- 
subevent method (i.e. squares in Figure 6) in the un- 
folding. This procedure could be employed to study 
the long-range correlations (known as the "ridge") re- 
cently observed in high multiplicity proton-proton [37] 
and proton- lead [38-40] collisions at the LHC. 

We appreciate valuable comments from R. Lacey. This 
research is supported by NSF under award number PHY- 
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